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Abstract This paper discusses a novel probabilistic approach for the design of robust 
model predictive control (MPC) laws for discrete-time linear systems affected by paramet- 
ric uncertainty and additive disturbances. The proposed technique is based on the iterated 
solution, at each step, of a finite-horizon optimal control problem (FHOCP) that takes into 
account a suitable number of randomly extracted scenarios of uncertainty and disturbances, 
followed by a specific command selection rule implemented in a receding horizon fashion. 
The scenario FHOCP is always convex, also when the uncertain parameters and distur- 
bance belong to non-convex sets, and irrespective of how the model uncertainty influences 
the system's matrices. Moreover, the computational complexity of the proposed approach 
does not depend on the uncertainty /disturbance dimensions, and scales quadratically with 
the control horizon. The main result in this paper is related to the analysis of the closed 
loop system under receding-horizon implementation of the scenario FHOCP, and essentially 
states that the devised control law guarantees constraint satisfaction at each step with some 
a-priori assigned probability p, while the system's state reaches the target set either asymp- 
totically, or in finite time with probability at least p. The proposed method may be a valid 
alternative when other existing techniques, either deterministic or stochastic, are not di- 
rectly usable due to excessive conservatism or to numerical intractability caused by lack of 
convexity of the robust or chance-constrained optimization problem. 



1 Introduction 

In Model Predictive Control (MPC), at each sampling time t, a plant's control input u t £ R m is computed 
by solving a constrained finite horizon optimal control problem (FHOCP), according to a receding horizon 
(RH) strategy, see, e.g., pQ. MPC has received an ever- increasing attention in the last decades, mainly 
due to the possibility of taking into account input and state constraints explicitly in the control design. 
The study of robust MPC approaches, able to guarantee stability and constraint satisfaction also in 
the presence of uncertainty/disturbances, is still a very active research area. For the case of linear time 
invariant (LTI) discrete time models, an extensive literature has been developed, considering the presence 
of either model uncertainty or external disturbances, see, e.g., [2]-[l2]- Most of the existing approaches 
are deterministic and aim to optimize a worst-case performance index, while enforcing constraints for 
all possible outcomes of the uncertainty [2]-[l] or disturbance [5]- [8]. These techniques guarantee that 
the designed control law is able to cope with the considered uncertainty. However, they rely on the 
assumption of convexity of the optimization problem, not only with respect to the control input, but also 
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with respect to the uncertain parameters and disturbances. Moreover, the computational complexity 
of deterministic approaches typically grows with the complexity of the model set. In a recent and 
active research direction, stochastic MPC techniques have also been studied, see, e.g., [9]- [12] and the 
references therein. Stochastic MPC techniques exploit some known statistical description of the uncertain 
parameters and/or of the disturbance (e.g., the probability distribution, or the first and second moments), 
yet they still employ deterministic algorithms and, in order to maintain tractability of the optimization 
problem, they typically assume that the system matrices are either perfectly known, or they have a 
particular structure that preserves convexity. 

We propose here a new randomized method for robust MPC design, which is able to deal with both model 
uncertainty and additive disturbances. Similar to stochastic MPC techniques, we exploit information on 
the statistics of the uncertain parameters and disturbances. However, we do not assume convexity or 
even connectedness of the model set or of the disturbance set. Still, the optimization problem in our 
approach is always convex, and the control law is able to robustly enforce constraints and trajectory 
convergence, with a probability higher than a user-defined value p. Furthermore, for a given value of 
p, the computational complexity of our approach is completely independent of the complexity of the 
model set. The key point enabling to achieve these features is a shift of paradigm, from a deterministic 
algorithm to a randomized one, i.e., an algorithm that relies on random choices in the course of its 
execution (see, e.g., [13]). Indeed, a key step in our main algorithm (Algorithm 14. lj) is the solution of a 
scenario FHOCP, in which we do not consider all possible outcomes of uncertainty and disturbances, but 
only a finite number M of randomly chosen instances of them, named the "scenarios." A randomized 
approach for MPC has been studied also in [141 ITS] , by using a Monte Carlo technique. However, Monte 
Carlo approaches may be very computationally demanding and can not handle in a straightforward way 
the presence of state constraints. Randomization has been used also in [16], in the context of chance- 
constrained MPC. However, in [16] there is no guideline on how to choose M in order to have the 
guarantee that the probability of success is at least p (which is instead one of the features of the present 
approach) and, moreover, the resulting optimization problem is a mixed-integer linear program. On the 
contrary, the approach proposed here, named MPCS (MPC via Scenario optimization), exploits relatively 
recent results in Random Convex Programming (RCP, see [17j-|20]) to provide an explicit link between 
M and p. Moreover, we introduce a slack variable, the "constraint violation," which renders the scenario 
FHOCP always feasible, and that can be used to monitor the extent of the (possible) violation of the 
involved constraints. Further, we show how scenario optimization can be embedded in a receding horizon 
scheme, in order to provide a feedback controller that gives probabilistic guarantees of robust stability 
and constraint satisfaction. The approach here proposed shall be particularly interesting in all those cases 
where the assumptions underpinning the existing deterministic or stochastic approaches for robust MPC 
are not met; for example, when the dependence of the system matrices on the uncertain parameters is 
not affinc. 



2 Problem formulation 

Consider the following uncertain, discrete time LTI model: 

x t+ i = A(6)xt + B(6)u t + Rj(0)it (1) 

where t £ Z is the discrete time variable, Xt £ K" is the system state, Ut £ K m is the control input, 
7t £ r C W 1 " 1 is an unmeasured disturbance vector, 9 £ Q C R 9 is the vector of uncertain parameters, 
and A(8), B(9), B 7 (6) are matrices of suitable dimensions. Let us consider the following assumptions: 

Assumption 1 (Uncertainty description) The sets T and £ = {A(6), B(9), B 1 (Q) : £ 0} are bounded. 
We assume "ft cind 6 to have stochastic nature, and we let Wg denote the probability measure on 0, and 
P 7 the probability measure on T. Variables and jt a re independent. Moreover, 7 = {70,71, • ■ •} is an 
independent identically distributed (i.i.d.) sequence and we let P2° denote the probability measure on this 
sequence. ■ 

Assumption 2 (Robust stabilizability) The pair A(9), B{9) is stabilizable for any 9 £ O. ■ 
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The control problem is to regulate the system state to a neighborhood of the origin, subject to (possibly 
uncertain) input and state constraints xt £ X(0), Ut £ U(0), Vi. The next assumption characterizes the 
constraint sets. 

Assumption 3 (State and input constraints) For any 9 € Q, the sets X(0) C R™ andV(9) C R m are con- 
vex; they contain the origin in their interiors and they are representable by: X(#) — {x £ R" : fx{x, 9) ^ 0} , 
U(0) ={»£ R m : fu(u,0) ^ 0}, where ^ denotes element-wise inequalities, each entry of the functions 
fx ■ R" x -> l r , /(/ : R m x — > R 9 is convex in x and u, respectively, and r, q are suitable integers. 

m 

The parameter 9 has been included in the constraints to account for practical applications where, 
for example, a convex function of the states (e.g., energy, load) has to be limited below some threshold, 
and some parameter in the function or the threshold itself depend on uncertain physical quantities (e.g., 
maximal energy, breaking load). Assumptions Q] and [3] are quite mild, since O may be unbounded and of 
any form, no assumption on S, T is made except for boundedness, no restrictions on how the parameter 
9 influences matrices A(6), B(9), B 7 (9) are imposed, as long as the system is stabilizable, and finally no 
assumption on the shape of the convex sets X(0), U(0) (e.g., polytopic, ellipsoidal, ...) for given 9 e is 
made. Mixed constraints of the form (x,u) £ X[/(0), where X;j(#) C 1" x R m is a convex set, are not 
considered here for simplicity, but they can be straightforwardly included in our problem settings. Due 
to the presence of the generally non-zero unmeasured disturbance 7 t , regulation of the system state to 
the equilibrium x = 0, u — can not be attained. Rather, we can require regulation to a neighborhood 
of the origin, described by a terminal set, which is robustly positively invariant under a terminal control 
law. 

Assumption 4 (Terminal set and terminal control law) A convex set X/ and a linear terminal control 
law u = KfX, K f £ R mxn ; exist for system (QJ), such that X f = {x : f Xf (x) ■< 0}; A(9)x + B{6)K f x + 
B 7 (0)7 £ X/, V6> £ 9, V 7 € T, Vz £ X /; finally f x (x,6) r< 0, fu{K f x,9) r< 0, V0 £ 6, Vz £ X /; w/iere 
/ X/ : R" -> R' has convex components, and I is a suitable integer. ■ 

A possible method for constructing a terminal set and a terminal control law satisfying Assumption |4] is 
to apply results in quadratic stability and rejection of bounded disturbances for uncertain LTI systems, 
see, e.g., [U] [52] and the references therein. Moreover, there is a number of contributions in the literature 
concerned with the computation of approximations of the (minimal) robust positively invariant terminal 
set X/, see e.g. [23l El] and the references therein. In the rest of this note, we parameterize the control 
input as: 

ut = KfX t + vt, (2) 

where Kf is the terminal control law of Assumption 2] (which is assumed to be known and given), and Vt 
is a control correction to be designed. Plugging filfy into ([J), we obtain the discrete-time model 

x t+1 = A cl {9)x t + B{9)v t + £ T (0)7 t , (3) 

with A c \(9) = A(9) + B(9)Kf, which will be the basis of our developments. 

3 The Scenario-based Finite-Horizon Optimal Control Problem 

Suppose that, at a given time instant t, the state xt of system (j3]) is observed. We consider the problem 
of determining a corrective control sequence on a horizon of N instants forward in time. To this end, 
we build a randomized finite- horizon optimal control problem (FHOCP), as described next. Let N be 
the chosen horizon length, and let Vju, j = 0, 1, . . . , N — 1, be the N predicted control corrections to be 
applied to ([3]), from t to t+N — 1, given the knowledge of the state at time t. From (j2J, the corresponding 
predicted control input sequence is Uju = KfXj\ t + Vj\ t , j — 0, 1, . . . , N — 1. By using model (j3]), we 
thus obtain the predicted values of the states as linear functions of the current state Xt, of the predicted 
(to-be-determined) control sequence Vt = [vZ t ■ ■ ■ J T £ R Nm , and of the disturbance sequence 7: 

x jlt = A j J6)x t + $,(0)Vt + T,(0) 7 , j = 1, . . . , N, (4) 
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where $j(6), Tj(0) are suitable functions of the model matrices, A c \(8), B{9) and B 1 {9). However, the 
predictions obtained via model ([3]) are uncertain, since they depend on 9 and on 7. In our approach, we 
deal with this issue by considering a discrete set of predicted state and input trajectories, obtained for 
a number M of randomly extracted scenarios of 9 and 7 at time t. More precisely, let us collect these 
random parameters in S = (9, 7), 5 G A = 6 x r°°. As a consequence of Assumption [TJ we have that S 
has a probability measure that we denote with P, which is the product measure of ¥g and the measure 
on 7: P = x P^°. Consider then M independent extractions £j , . . . , S[ M of 6, constituting the 

scenarios, where each scenario has the probability distribution P, and let uj t = (S^ 1 , . . . , §\ M ^ ) denote 
the "multisample" of scenario extractions at time t. The probability distribution on uj t is given by P M . 
Based on the random scenarios, we obtain M different state and input predictions from namely, for 
i = l,...,M, 

x o\t = Xt 

x% = 4((?W)x t + * i ((9W)V t + T i (^) 7 «, 
j = l,...,N, 



u<jj = K f xfl + v jlt , j = Q,...,N-l, 



where (0 t W , 7t W ) = 6^. Let us now introduce the following cost function 

f N-l N-l 

i=l,...,M 



J(x t ,u t ;V t ) = 5= max f ( £ d{x%X f ) + ]T vj lt Av jlt ] , (6) 

j'=o 3=0 



where d{x, X/) is the distance between x and the terminal set Xj , computed in some norm || • || : d(x, X/) = 
min \\x — y\\, and A = A T >- is a weighting matrix chosen by the control designer. In the following, 

yeXf 

with a slight abuse of notation, we indicate the state and input constraint sets as X(<5), U((5), respectively, 
and the related convex functions in Assumption [3] as fx{x, 6), fu(u, d). Moreover, we transform the hard 
constraints of Assumption [3] into soft ones, by introducing a slack variable q t G M, qt > 0. Then, the 
scenario-based FHOCP is a random convex program defined as follows: 

P(x t ,U t ) : min z t + aq t (7a) 

V t ,z t ,qt 

subject to 

J(x t ,u>t;V t ) < z t (7b) 

fx{x®,8®)-lqtl0; 3 = 1, ■ • ■ , N - 1, i = 1, . . . , M (7c) 

Muf^S^) - l ft d 0; 3 = 0, • ■ • , JV - 1, * = 1, • • • , M (7d) 

fx f (x?l mt )-lq t ±0; i = l,...,M. (7c) 

qt > (7f) 
In (Tfal) . the weighting scalar a > is chosen by the control designer, and 1 denotes a column vector of 
appropriate length, containing all ones. We denote with Vt(xt,0Jt) = { w o|t' * ■ • > v N-i\t }> ^(xti^t) and 
q^(xt,uJt) an optimal solution to problem (J7J- 

Remark 3.1 (Worst-case cost and constraint violation) Due to the presence of constraint (|7bl) . i/ie va/we 
z t * is an upper bound of the worst case cost with respect to all the M extracted scenarios. We thus refer 
to z t * as the "worst-case cost." Moreover, we note that the use of the soft constraints (|7c[) - ([7e] ) imply that 
problem (JT]) is always feasible. In particular, by using a sufficiently high value of a (e.g. 10 4 times higher 
than the typical value of z%), the optimal value of q^ turns out to be negligible whenever the problem with 
hard constraints (i.e., with q t set a priori to zero) is feasible. Contrary, when the problem with hard 
constraints is not feasible, the variable q^ provides an indication on "how much" some of the constraints 
are violated. For this reason, we refer to q% as the "constraint violation" level. Finally, we note that there 
is no constraint violation in (|7b[) . i.e. z t * is always greater than all the cost functions corresponding to 
the sampled scenarios, and in particular it is always an upper bound of the distance between the state Xt 
and the terminal set X/ (see ). This feature is important for our convergence result in Section^ ■ 
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Remark 3.2 (Choice of cost function and input parameterization) Prediction of the state trajectories 
in a closed loop fashion is quite common in the context of robust MFC, see e.g. |3 251. In particular, 
we adopt here the input parameterization |Q|), and we optimize over the control corrections Vj\ t , j 

0. Y 1 . i.e. over N m decision variables. Moreover, we chose as stage cost the distance between the 

state and the terminal set, plus a quadratic penalty on the control correction. Indeed, these choices of 
control parameterization and cost function are not meant to be the sole possibility, neither the optimal, 
for the proposed approach. Generalization to other kinds of input parameterization (e.g. disturbance- 
feedback {8, 251) and cost function (like a standard quadratic stage cost) can be done with some technical 
modifications in the proofs of the results reported in this note. ■ 

The optimization problem V(xt,ujt) can be rewritten in a more compact standard form. By collecting 
the optimization variables {Vt,zt, qt) in vector s t g R mN+2 , the cost can be expressed as z t + aqt — c T s t , 
where c = [0, . . . , 0, 1, a] T . Moreover, it can be noted that, for any fixed value of d t , due to linearity of 
([S]), the constraints (|7b|) - (|7fj) are convex in the decision variable st and in the state Xt- Finally, these 
constraints can be formally expressed compactly as h(st,x t ,S^) < 0, for all i — i,...,M, where h : 
R mN+2 x 1™ x A ^ M is defined as h(s t ,x t , S t (i) ) = max ( max \.fx(x% 5 t w ) - lq t , fu(u% 5®) - lq t \ , 

I j=0,...,N—l I ^ J 

fx f (»t+w|t) — — J{ x t-> u t'i Vt) — zt\ Notice that h(st,xt, 5^) is convex in both s t and Xt, since it is 
the point-wise maximum of convex functions. The scenario FHOCP can hence be rewritten as 

V(x t ,io t ) : min c T s t (8) 

St 

subject to: h(st,xt, 5^) < 0, i = 1, . . . , M. 

We denote with s^(xt,uJt) — (V ( * , z* t , q£) an optimal solution of V(xt,wt)- Notice that, due to the way 
it has been defined, problem V(x t ,uit) is always feasible. We further assume that this problem always 
attains a unique optimal solution. 

3.1 Properties of the scenario FHOCP 

We now consider the following problem: suppose that, given the state Xt, we solve problem V{xt,uJt)- 
Then, we ask what is the probability that the computed optimal control sequence V^{xt,ujt) is able to 
satisfy all state and input constraints over the chosen horizon, and to drive the state trajectory to the 
terminal set at the end of the horizon, within the computed optimal constraint violation g t * . Formally, this 
is the probability (with respect to S) with which h(s^ , xt, 5) < 0, where we notice that h is now evaluated 
at the optimal scenario solution s^, and the state and input trajectories that enter the definition of h 
are the "actual," uncertain, ones, obtained from model ([3]) at a random 5 — (0,j). So, we define the 
reliability R of the scenario-FHOCP as 

R = ¥{6 : h(s$,x t ,5) < 0}. 

Notice that R £ [0, 1] is itself a random variable, since it depends on s^, which in turn depends on 
the random multi-extraction of the scenarios uj t , hence R = R(uj t ). Indeed, for some extractions ujt the 
reliability can be good (close to one), and for other extractions it can be bad. It is therefore critical to 
assess the a-priori likelihood of these two situations, that is to precisely quantify bounds on the probability 
of the "bad" event where {R < p}, being p some a-priori assigned level or desired reliability. To this 
purpose, we exploit the fact that problem V(xt,uJt) belongs to the class of so-called Random Convex 
Programs (RCP) (see, e.g., [T7]-[20,) and, in particular, the result in Theorem 1 of [T5], concerned with 
feasible random convex programs, applies to our context. The following key result directly follows from 
Theorem 1 of [T^], see also Theorem 3.3 in [20] . 

Theorem 3.1 Let d — mN + 2 be the number of decision variables in problem V(xt,UJt), let p £ (0, 1) 
be a given desired reliability level, let (3 £ (0, 1) be a given small probability level (say, f3 — lO -9 ^, and let 
M be an integer such that 

$(jp,d,M)<0, (9) 
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d-1 /i.A 

with $(p,d,M) = E ■ (l-p) j p M - j . Then, it holds that 

3=0 V 3 J 



{u t : R(u t ) >p}>l-0. 



(10) 



Remark 3.3 (Number of scenarios and "certainty equivalence") The practical importance of the result 
in Theorem \3.1\ stems from the fact that the number M of scenarios necessary to fulfill condition 
grows mildly with the inverse of (3. More precisely, Corollary 5.1 in [20] states that condition {IJ) is 
implied by M > (ln^ 1 + d) , thus M grows at most logarithmically with (3 1 ( tighter values of M 
for given (3 and p can be obtained by inverting numerically ©J. This means in turn that the parameter 
(3 may be fixed by the designer to a very low level, say (3 = 10 -9 , or even (3 — 10~ 12 , and still the number 
M of scenarios necessary to guarantee \10\) remains manageable. The parameter f3 hence measures the 
probability of the unfortunate event in which the optimal solution has reliability smaller than the desired 
level p. If the likelihood of such an event is bounded a priori by an extremely low value, such as j3 = 10~ 12 ; 
then we may safely say that, to all practical engineering purposes, the event {R(co t ) > p} is the "certain" 
event. In other words, the possibility that {R(io t ) > p} is not satisfied by the scenario problem is so 
remote that, before having any concern about it, the designer should better verify the validity of many 
other assumptions and approximations in the model. Such a "certainty equivalence" principle, which will 
be adopted henceforth in this paper, essentially eliminates from consideration the "outer" probability level 
in M0\) . and states with practical certainty (the expression "with practical certainty" shall be used in the 
rest of this note as a synonym of "with probability larger than 1 — f3, " where j3 > is some extremely small 
value) that {R(uj t ) > p} holds. This simplifies greatly the practical application of scenario techniques, and 
makes the whole approach more clear and understandable by both theoreticians and control practitioners. 
■ 

The properties of the scenario FHOCP are resumed in the following proposition. 

Proposition 3.1 (Finite horizon robustness) Given the state Xt of system (0) at time t, consider the sce- 
nario problem V(%t,(jJt) as an instrument to derive a finite-horizon control sequence V t * = {v^ t , ■ ■ ■ , v N-i\t} 
to be applied to the system (Q) at the subsequent instants t, t+1, . . . ,t + N — 1. Let the number M of 
scenarios in problem V(xt,iOt) be chosen so to satisfy (G|) for given reliability level p 6 (0, 1) and very 
small (3 6 (0, 1). Then, with practical certainty it holds that the computed control sequence: 

a) steers the state of system to the terminal set Xy in N steps with probability at least p and constraint 
violation g t *, i.e.: P{<5 : fx f {%t+N ,o~) — l?t ^ 0} > p; 

b) Satisfies all state constraints with probability at least p and constraint violation q t * , i.e.: ¥{5 : 
f x (x t+j ,6)-lqt±0,Vje[l,N]}>p 

c) Satisfies all input constraints with probability at least p and constraint violation qf , i.e.: P{<5 : 
f u (u* tHl S)~lq^0,y J e[0,N-l}}>p. ' ■ 

The proof of this result follows immediately from Theorem 13.11 eq. (fT0|) states that, with practical 
certainty, the optimal solution sjf of the scenario problem satisfies h(s^ , Xt, 5) < with probability at 
least p, which indeed implies that points a)-c) in the corollary hold. 

Remark 3.4 (Relationship with deterministic approaches) In a deterministic approach to robust MPC, 
a problem similar to ([7]) has to be solved for all possible values of S S A. When the problem is convex with 
respect to S (which happens, for instance, when the uncertain matrices and/or the additive disturbance 
belong to polytopes f^\^j), deterministically robust approaches are indeed well-established and shall be 
preferred to the scenario approach, especially if deterministic robustness is critical in the considered 
application. In all other cases, deterministic approaches are generally intractable, unless the problem 
is manipulated so to satisfy convexity assumptions, at the cost of higher conservativeness and reduced 
feasibility. In these situations, the scenario approach proposed here is a viable alternative to deterministic 
techniques, since it is always convex and it can be efficiently solved also with a large number of samples, 
while still giving probabilistic guarantees on the robustness of the solution, as it will be shown in the 
example section. Moreover, we note that, if a sufficiently high weight a is used in the objective (see 
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Remark \3.1\) , then the scenario problem typically returns a negligible optimal constraint violation g t * , 
whenever the problem would be feasible with qt set to zero. Furthermore, since the scenario FHOCP has 
only a subset of the constraints of a corresponding deterministically robust FHOCP, the violation level 
ql of the scenario problem will always be lower than the violation level of the deterministic version of 
the same problem ( and this fact holds independently of whether the deterministic problem can be solved 
numerically or not). In any case, the value of g£ gives an indication on the extent of the violation of 
the involved constraints, which can be used, e.g., to implement supervisory control strategies and recovery 
actions. ■ 

The remaining part of this note is devoted to analyzing what happens when a scenario FHOCP is solved 
repeatedly in time and used to control the plant in a receding-horizon fashion. In a receding-horizon 
approach, which is the key feature of MPC, only the first control correction in the optimal sequence V t * is 
applied at time t, and then the FHOCP is solved again at time t + 1, by exploiting the knowledge of the 
state Xt+i, etc. In the next section, we propose a technique for incorporating the scenario FHOCP into a 
suitable receding-horizon scheme, and we derive probabilistic guarantees of asymptotic convergence and 
constraint satisfaction for the resulting closed-loop system. 



4 MPC scheme based on Scenario optimization 

We here introduce a receding-horizon implementation of a control algorithm based on the scenario 
FHOCP, as described next. The notation is set as follows: "*" variables, such as z t *, q% , V t * = {v% t , • ■ • , v N~i\t}> 
denote the optimal solution of the scenario optimization problem V{xt,uit) at time t, given xt\ "~" vari- 
ables, z t ,qt,Vt, denote, respectively, two scalar values and a sequence of N vectors of dimension m, as 
defined in the algorithm below; finally plain variables, Zt,qt,Vt, denote the running values of the variables 
z, q and of the sequence V = {«o|t> • • • > v N-i\t} in the algorithm. The first entry in Vt, namely Vo\ t , is the 
actual control correction that is applied to the system ([3]) at time t. The subsequence composed by the 
last N — 1 elements of Vt is denoted with ui : jv-i|t- We are now in position to describe the algorithm for 
MPC based on Scenario optimization (MPCS). 

Algorithm 4.1 (MPCS algorithm) 
(Initialization) Choose a desired reliability level p £ (0,1) and "certainty equivalence" level /3 £ (0,1) 
(say, (3 — 10~ 9 , or j3 = 10~ 12 J. Let M be an integer satisfying fP|). Choose e £ (0,1] (see Remark \4-.1\ 
below for the meaning of e and for guidelines on its choice). Given an initial state xq, extract ujq according 
to P , solve problem Vm(xo,ujo) and obtain the optimal control sequence Vq = {Wqi , v*< . . . , w^-_ 1 | }, 
and the optimal objective Zq and constraint violation <7q. Set Zq = z^,qo = q$, Vo = Vq, and apply to the 
system the control action uq — KfXQ + v \ . 

1) Lett := t + 1, observe x t , and set V t = {«i|t_i, . . . , Wjv-i|t-i>0} = { v v.N-i\t-i, 0}, 
z t = max(0,zt_i - d(aj t _i,X/)) , q t = q t -i; 

2) Extract the multi-sample Lu t according to P M , and solve problem VM{xt,uJt). Let (V t *, z*, ql ) be the 
obtained optimal solution. 

3) Evaluate the following collectively exhaustive and mutually exclusive cases: 

3. a) If z*[. > (zt_i — ed(xt-i, X/)) and St < d(xt,X/), then set Vt — Vt; z t — 0; qt = qt,' 
3.b) If Z? > (z t -i - ed(x t -i,Xf)) and z t > d(x t ,Xf), then set V t = V t ; z t = z t ] qt = qt,; 
3.c) If Zt < (z t -i - ed(x t -i,Xf)), then set V t = V t *; z t = z^; q t = q%; 

4) Apply the control input Ut = KfXt +fo|t> then go to 1). ■ 

Remark 4.1 The inequality z\ < {zt-\ — sd(xt-i, X/)), checked at step 3) of the MPCS algorithm, can 
be interpreted as a verification of a required minimum improvement, in terms of worst-case cost, achieved 
by the newly computed optimal solution (V t *,^ t *,gj) of the scenario problem at time step t, with respect 
to the previous step. The user-defined parameter e £ (0, 1] influences such a requirement: the closer the 
value of e is set to 0, the more likely it is that case z\ < {zt— l — sd[xt— i, X/)) is met, so that the MPCS 
algorithm relies, at each time step, on the newly computed optimal solution. Vice-versa, the closer is 
the value of e to 1, the more likely it is that the complementary condition z\ > {z t -i — ed(xt-i, Xy)) is 
detected, so that the MPCS algorithm employs the previously computed solution. ■ 
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The next results is concerned with the guaranteed properties, in terms of constraint satisfaction and 
convergence to the terminal set, of the closed loop system obtained by applying Algorithm 14. II 

Theorem 4.1 (Properties of Scenario MPC) Let Assumptions]]]^ be satisfied and letp G (0, 1) be a cho- 
sen reliability level. Let Vqu, t = 0, 1, . . . denote the sequence of control actions produced by Algorithm ^. 1\ 
and consider the closed loop system obtained by applying to (QJ) the control law Ut — KfXt + «o|t- Then: 

(a) With practical certainty, at all time steps t = 0,1,..., the probability that the state and input 
constraints are satisfied with constraint violation qt is at least p, that is P{<5 : fx(xt+i,S) — l?t dt 
n fu(u t ,S)-lq t ^0}>p, t = 0,l,... 

(b) Algorithm ^. 1\ either: (i) makes the state trajectory converge asymptotically to the terminal set, i.e. 
lim d(xt,Xf) = 0, or (ii) there exists a finite time t* such that, with practical certainty, the control se- 

t— »OC 

quence {t>o|t* i u o|t*+i> ■ • ■ v o\f+N-i} drives the state of the closed-loop system to the terminal set at time 
t* + N — 1, with probability at least p and constraint violation q t * ■ I 

Proof 4.1 Preliminaries. Notice first that, for any t > 0, if x t 6 X^ then the optimal solution to 
problem V(xt,uit) is (V*,z*,q^) = (0,0,0), since the terminal control law (i.e., with Vt = 0) is able to 
keep the predicted state trajectory in the terminal set while satisfying all constraints. Also, if z t * = is 
the optimal objective of problem V(xt,Uit)> then xt G Xy, since z t * is an upper bound o/d(xt,X/) (see also 
Remark 3.1), therefore z 4 * = x t G X/. Let then xt £" X/. 

Proof of statement (a). At time t — 0, Provosition \3.1\ guarantees with practical certainty that the 
first control correction satisfies the constraints on uq and x\ with probability no less than p and constraint 
violation q t = q*. At any generic time step t > 1, the variables (Vt,Zt,qt) are computed. Then, two 
cases may occur. If z t * < (zt-\ — ed(xt—i, X/)), then case 3.c) is detected, and the first element v^ t 
of the optimal sequence V t * is applied to the system. Being this sequence the solution of a scenario 
optimization problem, with practical certainty the probability of satisfying state and input constraints is 
no less than p, with constraint violation q t — q% . If, on the other hand, z t * > {z t -\ — ed(xt-i, X/)), 
then we are either in case 3. a) or S.b), and in both cases the element v^ t _ k , for some k G [1,N — 1], is 
applied to the system. Being this value part of the solution sequence V^_ k , with corresponding constraint 
violation q^_ k , again the probability of satisfying state and input constraints is no less than p, with 
constraint violation qt — qt — a t-k- Thus, in any case, with practical certainty, at each time step the 
MPCS algorithm guarantees satisfaction of state and input constraints with probability no less than p and 
constraint violation qt- 

Proof of statement (b). Each run of Algorithm J^.l may have one of two possible behaviors, depending 
on whether or not there exists a finite time t > such that z* > (z t ~\ — sd(xt~i,^f)) and zt < d(xt,'Kf), 
that is, whether or not the situation in step 3. a is ever satisfied. We then name A the situation when 
condition in step 3. a is met at some finite t > 0, and A the complementary situation when this condition 
is not satisfied at any finite time, that is when z t * < (zt-i — ed(^_i,X/)) or z t > d(^t,X^) holds for all 
t>0. 

I. Let us first consider the situation of case A. Consider a generic time t. At step 3) of the MPCS 
algorithm, if z 4 * > (z t -i — ed(xt-i, X^)), then, since it is assumed that we are in situation A, it must 
hold that Zt > d(xt,^f), thus case S.b) occurs, and the values Vt — Vt and z t — z t are set. Now, recalling 
that Zt — max(0,2: t _i — d(xt-i, X^)), two cases may occur: either z t — or z t — Zt-i — d(xt-i,^f) > 0. 
If z t — 0, we have = z t > d(s(,Xf), i.e. d(xt,X/) = 0, which would imply that the terminal set has 
been reached. Otherwise, if z t = z t ~\ — d(x t _i,Xy) > 0, then we have: 

z t = Zt>d(x t ,Xf)>0, (11) 

and z t - z t -\ =z t ~ z t -i = z t -i - d(x t -i,Xf) - z t -i = —d(x t -i, X/). Thus, 

zt-zt-i < -ed(xt-!,X f ), Va: t _i &X f , (12) 

and Zt — Zt-i — Xt-i G X/. (13) 

On the other hand, if at step 3) of the MPCS algorithm it happens that z t * < (z t _i — ed(xt-\, X/)), then 
case S.c) occurs, and the optimal values Vt and z t * are retained, i.e. z t — z^ , Vt = V t *. In this case, it is 
straightforward to note that equations ill\) - ![T3\) still hold true. The same reasoning can be repeated for 
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any time step, as long as the case z t * < (zt-i — ed(xt—i, X/)) or z t > d(xt,Xf) holds true as assumed, 
so that we can conclude that the variable Zt enjoys the following properties: 

z t > d(x t ,X f ) > 0, V< > 

2* = o x t ex f ( u] 

z t+1 - z t < -sd(x t ,X f ), Vx t <?X f , Vi>0 1 ' 

zt+i - z t = i t e X/. 

Properties fi^| j are sufficient to prove convergence of the state x t to the set Xf: 

< lim X/) < lim z f = 0, =>■ lim <i(a:t, X/) = 0. 

t— )-oo i— >oo t— >oo 

Therefore, we obtain that in case A the MPCS algorithm guarantees that lim c?(x t ,X*) = 0. 

77. Lei ms nerrf analyze what happens in case A. Let t > be the time instant at which the case 
zl > (zt-i — ed(xt-i,Xf)) and St < d(xt,Xf) is met for the first time, and let t* < i be the last 
time at which case z* < (zt—i — ed(xt—i,Xf)) was satisfied, that is the last time previous to i when an 
optimal command sequence was retained, together with its constraint violation , according to case S.c) 
of Algorithm J^.l; let I = i — t* > 1. According to step 3. a) of the MPCS algorithm, we set 

Vt = % zt = 0, q- t = qt- (15) 

Thus, at step 4) of the algorithm, the control move uj — KfXt + v ^ is applied to the system at time 
i, where v ^ — v*^ t , , i.e., v ^ is the optimal correction predicted for time t* + £ — i, computed at time 

t* . At time step t = t+ 1. the state variable xi + i is observed and [Vj+i, <7t+i) o,re computed as 
zf+i = max(0, zi - d{xi, X/)), q i+1 = qi, Vf+i = {n-.N-i\i, 0} = {«| +1 | t , , ^ +2 |t* > ■ ■ ■ > v *N-i\f , 0, . . . , 0}. 
Since H5\) holds, it must be % +1 = 0. Then, z|, 1} and V| +1 are computed at step 2), and we 

notice that, by definition, z^ +1 > 0. Therefore, at step 3) of the algorithm either (i) case 3. a) {z^ +l > 
(zi — ed{xi,Xf)) and zi+i < ^(^t+ii Xf)} is detected again, or (ii) one of cases 3. b) or S.c) are detected, 
which would imply, respectively, = % +1 > d(xt+i,Xf), or < d(xt+i,Xf) < z? +1 = % +1 = 0. Hence 
(in either case) xt+i £ X/, so that convergence to the terminal set would be achieved. Consider then 
case (i): the values Vj +1 = Vj + i, = and qt+i = qi ar £ set in the algorithm, and the control move 
ui + i = KfXt + i+Vp +1 ^ t , is applied to the system. Now, the same circumstances actually reproduce for all 
time steps t = t + k, k > 0, so the algorithm is such that the optimal input sequence V** , computed at time 
t* by solving a scenario FHOCP, is the one actually next applied to the system, and the related constraint 
violation g ( * is retained for all t > t* . Thus, in case A, there exists a finite time t* such that the sequence 
V 4 * is applied to the system for all subsequent instants t = t* + k, k = 0, . . . , N — 1. Now, the sequence 
V t * is the result of the solution of the scenario- FHO CP V(xt- , Wt*)> and Provosition \3.1\ states that, with 
practical certainty, we have R(ujf) > p, where R is the reliability defined in Section III- A of the paper, 
which means that ¥{S : h(st* ,x t * ,S) < 0} > p. Therefore, in the situation A, there exists a finite time 
t* at which an optimal control sequence is computed by solving a scenario-FHOPC and next applied to 
the actual system for the subsequent N time instants: we can hence claim with practical certainty this 
sequence will satisfy the problem constraints and reach the terminal set within the time window from t* 
to t* + N , with probability at least p and constraint violation g t *. 



5 Numerical example 



We consider the system (Q} with: A(8) 



1 + 8 
O.lsinf 



0.3 arctan(#5) 
1 



l + 6> 3 



1 + 6! ,B{0) 

7 4 j 1 + 6*2 

Each of the random parameters 9±, O2, #3 is uniformly distributed in the interval [-0.1,0.1], while #4, #5 
are distributed according to Gaussian distributions with zero mean and unit variance. Moreover, the 
disturbance 74 g M 2 is computed as follows: 



It = 



{[ 

I [C 



ni, min (772, g ( 100(3^1+0.05) - 
0.05 cos (773) , ?7 5 ] T if 770 < 0.5 



0.05 



))]' 



if 770 > 0.5 



9 



where 



r/ 5 = max ( min I 7/4 sin 



J), 0.05| sin (ijs)|), -0.05| sin (773) |) 



and 770 € [0,1], rji € [0,0.05], 772 £ [0,0.05], 773 £ [§7T, §7r], 774 G [-0.05,0.05] are uniformly distributed 
random variables. We also consider the following constraints on the input and state variables: X(0) = 
{iel 2 : Halloo <x{0)} , U(0) = {u £ R : |u| < , where u(0) = 5/(1 + 6> 6 sin(0 7 )), z(0) = [10/(1- 
6 sin(0 7 )), 10/(1 + 6» 6 cos(6» 7 ))] T , 6> 6 G [-0.05,0.05] is uniformly distributed, finally 9 7 is distributed 
according to a Gaussian distribution with zero mean and unit variance. Note that the system matrices 
and constraints are a nonlinear functions of the uncertain parameters, and the disturbance belongs to a 
non-convex, disconnected set. Hence, no existing technique for robust MPC can be directly applied in 
this example. By using the results in [3TJ[32] (see Section [5]), we computed the terminal control law Kf 
and terminal set X/ satisfying Assumptions) K f = [-0.4686, -1.4221], X/ = {x G R 2 : x T Q f x < 1}, 



Qf = 0724 1724 ' We desi S ned the Mp CS law with N = 10 and A = 1. We set f3 = 10 J , and 
considered different values of p. In order to estimate the probability with which the MPCS control law 
and the finite horizon sequence Vq satisfy the constraints and drive the state to the terminal set, and to 
compare it with the bound p, we carried out Aerials = 100, 000 Monte Carlo simulations, starting from the 
same initial state xo = [5, 2.75] T . We note that this initial condition is not feasible for the deterministic 
counterpart of the scenario problem, hence for some extractions of loq the constraint violation Qq is not 
negligible. In each one of these simulations, Algorithm 14. ll has been applied and the probability of success 

p has been estimated as p = — trials failures ^ wnere _/Vf a ii ures is the number of simulations in which some 

Atrials 

of the constraints were not satisfied. Moreover, the finite horizon solution was also tested, to check the 
result of Proposition 13.11 In particular, for the finite horizon sequence Vq, the considered constraints 
were the state and input constraints X(0), U(0), and the terminal set constraint xn G X/, while for 
the receding horizon implementation the latter constraint was replaced with xjv +1 q G X/, in order to 
approximately take into account the asymptotic stability result of Theorem 14.11 (b)-(i). We evaluated 
these constraints as hard constraints, i.e. with zero constraint violation. The values of p for the finite 
horizon simulations and for the receding horizon ones are indicated as p FH and p RH , respectively. The 



Table 1: Numerical example. Estimates p of the probability of success for different values of p, with 
/3 = 10~ 9 . 



M[p) 


22 (0.05) 




42 (0.30) 






P FH = 0.885, p Hti 


= 0.921 


p h ti = 0.901, j5 HH = 


0.943 


M(p) 


95 (0.60) 




890 (0.95) 






p b ti = 0.923, p Hti 


= 0.963 


p b ti = 0.993, p Kti = 


0.999 



obtained results are reported in TableQ]for values of p = 0.05, 0.3, 0.6, 0.95. The corresponding values of 
M are also reported in the Table. It can be noted that in all cases the values of p FH , p RH are higher than 
the corresponding p, in accordance with the theoretical results of Sections [SJH] Moreover, the values of 
p FH are lower than those of p RH . The estimated probabilities of success p FH , j5 RH are quite good already 
with low values of M: in the finite horizon case, the reason of such a good result is mainly the presence of 
the terminal control law Kf, which has already some degree of robustness, while in the receding horizon 
case, a higher robustness derives from the iterative re-optimization of the corrective control sequence V t *. 
It is worth to notice that the value of M does not depend on the dimension of the state variable or of 
the uncertainty/disturbance variables; it only depends on the chosen values of p, /3 and on the number of 
decision variables in the scenario FHOCP, i.e., the number m of inputs multiplied by the control horizon 
N, plus the slack variables z and q. However, the number of constraints embedded in h(s,x,S) depends 
linearly on n, m and N, so that the growth of the overall number of constraints in the scenario problem, 
for fixed values of p and /3, is ~ (n • m 2 ■ N 2 ), i.e. quadratic with respect to the control horizon. 
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